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In this paper we introduce a method of calculating the local temperature and chemical potential 
inside a mesoscopic device out of equilibrium. We show how to check the conditions of local thermal 
equilibrium as the whole system is out of equilibrium. Especially we study the onsite chemical 
potentials inside a chain coupled to two reservoirs at a finite voltage bias. In the presence of 
disorder we observe a large fluctuation in onsite chemical potentials, which can be suppressed by 
the electron-electron interaction. By taking average with respect to the configurations of disorder, 
we recover the classical picture where the voltage drops monotonously through the resistance wire. 
We prove the existence of local intensive variables in a mesoscopic device which is in equilibrium or 
not far from equilibrium. 
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I. INTRODUCTION 

In recent years the thermodynamic properties of a 
quantum system driven out of equihbrium have attracted 
more and more attention [1-7]. Especially, in meso- 
scopic transport a standard method of studying nonequi- 
hbrium steady states (NESS) has been established b 
using the evolution approach in Keldysh formalism 
[To| . In mesoscopic transport we study small quantum 
systems with a few degrees of freedom which must be de- 
scribed by quantum mechanics. The small system is con- 
nected to several infinite reservoirs which are in thermal 
equilibrium with different temperatures and chemical po- 
tentials. The nonequilibrium steady state is approached 
by an evolution starting from the supposed initial con- 
dition. The evolution approach has been widely used 
to study the transport through quantum point contacts, 
quantum dots, single molecules and carbon nanotubes. 
Good agreements between theoretical predictions and ex- 
perimental results have been obtained. 

Now a fundamental question arises as to whether the 
concepts in equilibrium thermodynamics apply to the 
mesoscopic system driven out of equilibrium. For ex- 
ample, temperature is an important concept in equilib- 
rium thermodynamics, where it is also called the inten- 
sive thermodynamic variable of a system. According to 
the zeroth law of thermodynamics, it is always possible 
to assign a temperature to a system in thermal equilib- 
rium. Even if the whole system is out of equilibrium, 
when the conditions of local equilibrium are fulfilled, we 
can divide the system into small cells which are in ther- 
mal equilibrium to a good approximation. In principle 
we are then able to define the local temperature in each 
cell. However, up to now the sole way of calculating tem- 
perature according to the density matrix of a system is by 
using the extremal principle in statistical mechanics. Ac- 
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cording to Glansdorff and Prigogine llj in systems away 
from equilibrium there is no general variational principle. 
Then temperature, arising from the Lagrangian multi- 
plier in statistical mechanics, loses its meaning. Even 
if many authors make their efforts in the local extremal 
principles [TT| - [23j which can be used for local solutions 
when the system is not far from equilibrium and dissi- 
pative processes dominate there, there is no consensus 
in if these principles can be applied to mesoscopic trans- 
port [Hi. 

The mesoscopic transport should have been a touch- 
stone of the principles in nonequilibrium statistical me- 
chanics. But it is not the case because of lack of defini- 
tions of local intensive thermodynamic variables, which 
play central roles in statistical mechanics. In evolution 
approach even if the currents can be easily calculated, 
up to now there is no way of calculating the temperature 
in a local region of the mesoscopic device. Because the 
current and the temperature are two different kinds of 
physical quantities. The current is calculated as the ex- 
pectation value of the quantum operator of current with 
respect to the density matrix of the system. While the 
temperature is intensive, which in statistical mechanics 
is introduced by hand as a Lagrangian multiplier before 
one writes down the density matrix. There is no "tem- 
perature operator" , although in fact temperature is also 
measurable like all the other physical quantities. Such 
a discrepancy between temperature and other physical 
quantities causes troubles in the theory of quantum trans- 
port. In fact there is no reason for ignoring the local in- 
tensive variables in mesoscopic device, since in principle 
we could continuously increase the device size without 
meeting a breakdown point of the physical laws behind, 
i.e. the Schrodinger equation, until the device is in the 
macroscopic region where the driving forces are the gra- 
dients of intensive variables rather than their differences 
in reservoirs. 

In this paper we get rid of this discrepancy by introduc- 
ing a method of checking if the local equilibrium condi- 
tion is fulfilled in arbitrary nonequilibrium steady states 
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and of calculating the local temperature and chemical 
potential. We apply our method in a disordered chain 
subject to a current, and study the onsite chemical po- 
tential by attaching an auxiliary site to it. We study 
the distribution of onsite chemical potentials at different 
disorder strength and reservoir voltage bias. In Sec. II 
we discuss how to reach the nonequilibrium steady state 
in an evolution approach by taking thermodynamic limit 
before taking t oo. In Sec. Ill, we discuss the condtion 
of local equilibrium and suggest a method of calculating 
the local intensive variables. In Sec. IV, we introduce the 
model of disordered ID chain and show the distribution 
of onsite chemical potentials inside it subject to finite 
voltage bias. Sec. V is a short summary. 



II. WHEN WILL A SYSTEM EVOLVE INTO A 
NONEQUILIBRIUM STEADY STATE 

When we say that a system evolves into a nonequilib- 
rium steady state, we mean that when time goes to infin- 
ity the system relaxes towards a stationary state in which 
some flows are nonzero, but there is no time variation. 
The NESS should be distinguished with the equilibrium 
state in which some nonzero flows exist in a ring struc- 
ture. The latter can exist in an isolated system with 
finite number of degrees of freedom, such as the meso- 
scopic ring in a magnetic field with a persistent current. 
But an isolated system with finite number of degrees of 
freedom will never evolve into NESS. Since if there exists 
a flow in an irreversible process, the flow carries some 
quantity (mass, energy or charge) from one part of the 
system to the other part, causing a continuous decreasing 
of the quantity in one part and increasing in the other 
part due to the conservation law. This contradicts the 
steady state argument. 

To avoid the paradox between irreversibility and sta- 
tionarity, we embed the system in an environment so that 
the composite isolated system is infinite. Then the total 
quantity in the whole system is infinite, so that any global 
conservation law is invalidated. At the same time if we 
see the local subsystem, the density of matter and the 
flow are both time invariant. 

In practice, to approach a nonequilibrium steady state 
we generally begin with a finite model and then increase 
the number of degrees of freedom into infinity. This pro- 
cess is called taking thermodynamic limit. One must no- 
tice that in most cases taking thermodynamic limit and 
taking t ^ oo are not exchangeable. To get the nonequi- 
librium steady state one must take thermodynamic limit 
before taking t — > oo. Next we give two examples, namely 
single impurity Anderson model and an infinite chain. 



A. Single impurity Anderson model 

Let us recall single impurity Anderson model which 
describes a quantum dot coupled to two leads (the left 



and right leads). This model is exactly solvable without 
considering electron-electron interaction and has been ex- 
tensively studied in condensed matter community. It is 
well known that a stationary current through the dot will 
be established if the two leads arc in different chemical 
potentials and the coupling between leads and the dot 
has been switched on in the infinite past. 

The crucial condition for a stationary nonzero current 
is that the number of levels in leads must be infinite. 
Otherwise the current perpetually oscillates with an av- 
erage value of zero. This is easy to see if we suppose there 
is only one level in each lead. Then the Hamiltonian is 
written as 

H = uj{c{cL + c^Cfl + d^d) + g{cld + c^j^d + h.c), (1) 

where cl, cr and d denote the annihilation operators in 
left lead, right lead and the dot respectively. The left lead 
is occupied by an electron at time t = 0. After switching 
on the coupling, we find the current to be 

lL{t) = ^\g\ sin (V2\g\t) + ^\g\ sin {2V2\g\t) , (2) 
which satisfies the conservation law 



Here Nl is the electron number in left lead. Whatever fi- 
nite number of levels there are in leads, the time-averaged 
current must be zero. Otherwise the current will even- 
tually empty one lead and overflow the other. In ther- 
modynamic limit, however, the term in righthand side of 
Eq.[3]is nonsense since is infinite. So the conservation 
law does not prohibit a nonzero current any more. 

One should notice that the stationary current in NESS 
is not carried by the eigenstates of the Hamiltonian. The 
invariant term of current operator in Anderson model is 
zero. In other words, the expecation value of II is zero 
with regard to arbitrary eigenstate, because the current 
operator is II ^ [H,Nl\- The fact that the eigenstates 
carry no current distinguishes NESS with equilibrium 
states carrying nonzero current. 

B. An infinite chain 

In single impurity Anderson model we employ the open 
boundary conditions, while the NESS can also be ap- 
proached in an infinite chain with periodic boundary con- 
ditions. 

The Hamiltonian of a chain of length 2N with periodic 
boundary conditions is written as 

2N-1 

H = -g + h.c.) ~ giclj^ci + h.c). (4) 

i=l 

In momentum basis it is diagonal as H = X]fc=o'^ ^kc\ck 
where efc = —2gcos{2Trk/2N). The operator at site j is 
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expressed in momentum basis as 

i2TTjk/2N 



E 

/c=0 



-Cfc. 



(5) 



Then we immediately have 

2N 



where Wj 
tor. 



2Af-l 
fc=0 



1 

i/v' 



is the propaga- 



Now we suppose that the sites from j — 1 to j — N 
are all occupied by electrons and sites from j = + 1 to 
J = 2N are empty at initial time t = Q. Then the hopping 
interaction in Eq.|3]is switched on. The local current and 
electron density at arbitrary time are worked out easily, 
given the current operator from site j to (j + 1) as 



(6) 



Obviously this definition satisfies the local conservation 

drii T 



law 



dt 



The current evaluates to 



N 



/,^,+i(t) = 2.9^ImW^*.,W^,+i,,v, 
and the electron density at site j evaluates to 



nj(t) 



N 



(7) 



(8) 



According to the behavior of scaled parameter g as 
iV — oo, the system will evolve into either a NESS with 
nonzero local current or a steady state with zero current 
everywhere. 

If g/N keeps a constant, the electron density nj{t) and 
the ratio Ij^j+i (t) / g oscillate around some values with 
a period which does not vary with N (see Fig. 1). As 
N ^ oo, nj{t) and Ij^j+i{t)/g both have well defined 
thermodynamic limits. For arbitrary j, the oscillation 
amplitudes of nj(t) and Ij^j+i{t) / g decay in course of 
time. When time goes to infinity the electron density 
approaches to 0.5 and the current to zero everywhere. 

If g keeps a constant, however, as N increasing the 
period of functions nj{t) and Ij^j+i{t) increases accord- 
ingly. For any finite N , the current eventually decays to- 
wards zero. But if we take thermodynamic limit firstly, 
the period of oscillations goes to infinity. The current at 
some sites will approach to a nonzero value when we take 
t — ^ CXI thereafter. 

Let us first study the currents at sites j = (1 + x)N as 
< a; < 1 is a constant. One should notice j will vary 
with A'^. In thermodynamic limit these sites are infinitely 
far away from the electron reservoir (the sites from 1 to 
N), so that the current 1 and the electron density 
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FIG. 1; The evolution of current from site A'^ to -I- 1 for 
different chain length when g/N is a constant. Here we set 
g = N/800. The y-axis is the ratio of current to g and x- 
axis the time. One should notice that the unit of time is not 
shown in the figure. But it is in fact 800N/g. The function 
lN^N+i{t)/g has converged when N is as large as 800. 
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FIG. 2: The evolution of current for different N when g is a 
constant. Here we set g — 1. With increasing A' the period of 
current will go to infinity. If one takes N ^ oo before taking 
t — >■ CX3, the steady limit of current will be around 0.636. 



Uj both approach to zero. On the other hand at the site 
j — {I + x)N as — 1 < X < 0, which is in the reservoir, 
the current approaches to zero and the electron density 
approaches to one. 

The site A^ -f n as n is an arbitrary integer is in 
the boundary between the reservoir and the vacuum. 
There the current approaches approximately to 0.636 (see 
Fig. 2), a nonzero value. At the same time the electron 
density approaches to 0.5. Here n can be arbitrarily large 
or small since in thermodynamic limit we always have 
\n\ N whatever n is. Then in the system consisting of 
sites N + n, which can be even macroscopic in size, there 
exists a well defined nonequilibrium steady state. 

As shown above, solving equations of motion and tak- 
ing thermodynamic limit before taking t — > oo supply a 
platform for studying the NESS. In next section we will 
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discuss how to define the thermodynamic intensive vari- 
ables, such as temperature and chemical potential, in the 
NESS got by this approach. The definition of intensive 
variabls is important both in experiments where they can 
be measured and in the theory of nonequilibrium ther- 
modynamics where they must be given before hand. 



III. HOW TO CALCULATE LOCAL 
TEMPERATURE AND CHEMICAL POTENTIAL 
IN A SYSTEM IN NONEQUILIBRIUM STEADY 
STATES 

Today quantum mechanics is believed to be the univer- 
sal theory describing the evolution of both microscopic 
and macroscopic systems. In quantum mechanics the 
observables are represented by linear operators acting 
on the Hilbert space. On the other hand the concepts 
in themodynamics, such as temperature T and chemi- 
cal potential /it, have been widely used to describe the 
macroscopic systems in everyday life. But unfortunately 
these intensive variables in thermodynamics cannot be 
directly related to quantum operators. Instead, in statis- 
tical mechanics T and /i are introduced by hand as the 
Lagrangian multipliers for the energy H and the parti- 
cle number N respectively. The lack of interpretations 
for T and /i at the microscopic level causes problems in 
nonequilibrium thermodynamics, where there is no gen- 
erally accepted variational principle and then the tem- 
perature and chemical potential (or more strictly the 
local temperature and local chemical potential) cannot 
be determined by the density matrix of the system. In 
other words there is no way of calculating T and /x in 
a NESS which is not the result of variational principle 
but the steady limit of the state evolving according to 
Schrodinger equation. We contribute to solve this prob- 
lem by providing a method of calculating the local tem- 
perature and chemical potential in a system according to 
its density matrix. 

To avoid ambiguity we would like to emphasize that we 
do not define local temperature and chemical potential in 
arbitrary NESS. In fact the definitions of T and /i make 
sense only in systems where the local thermodynamic 
equilibrium conditions are fulfilled. In other words, if we 
divide the system into small cells each cell must look like 
being in an equilibrium state, even if the whole system 
undergoes an irreversible process. In systems violating 
the local equilibrium arguments one would get arbitrary 
T and /i by employing different measurement schemes 
which give unique result in an equilibrium system. As 
will be shown, an approach examining if the local equi- 
librium conditions are fulfilled naturally arises from our 
method. 

Let us recall the process of using a mercury thermome- 
ter to measure the temperature of an object. We attach 
the thermometer closely to the target and wait for long 
enough until the temperatures of the target and the ther- 
mometer are the same. If the target is much larger than 



the thermometer, its temperature will not be changed 
during the equilibration process. Then we read out the 
temperature according to the volume of mercury. This 
procedure suggests a general scheme of measuring local T 
and /i of a system either in equilibrium or out of equilib- 
rium by attaching an auxiliary apparatus to the region we 
want to measure. After waiting for a period longer than 
the equilibration time, the auxiliary apparatus will equi- 
librate if the conditions of local equilibrium are satisfied 
even if the whole system may be still out of equilibrium. 
Then T and /x are read out in the auxiliary apparatus ac- 
cording to their orthodox definitions in equilibrium. The 
auxiliary apparatus (the "thermometer") should satisfy 
the following conditions: it must be very small so that it 
almost has no effect on the original state of the system; 
and there is a simple relation between its temperature 
and some measurable property. Obviously a single site is 
a suitable candidate of such auxiliary apparatus, since it 
is the smallest apparatus and its occupation number for 
fermions can be expressed as 

= e^(^.-p)+l ' 

where e^, /3 and /i are the energy level, the inverse of 
temperature and the chemical potential respectively. 

Our method bases on an ambitious assumption that a 
single site will equilibrate after attached to a large sys- 
tem. This is true only under special conditions, i.e., the 
coupling between the site and the system is infinitesimal. 
Consider an impurity site coupled to a Fermi sea with 
chemical potential ^ and temperature T . This is just 
the single impurity Anderson model without interaction. 
We switch on the coupling at time t = Q. Given the en- 
ergy level td and the level width T at the impurity site, 
in steady limit when time goes to infinity we find the 
occupation number to be 

'^^ = ^y_'V(^-) + i(.-..)^ + r- ^10) 

Only in weak coupling limit F — > wc have 
7(7377p+r^ 5{e — e^) and then recover Eq. [HI We 
rewrite Eq. [9] as 

In (^i- - 1^ = Pea - Pfi. (11) 

Then /3 and /i can be easily determined by plotting 
hi{l/nd — 1) with respect to e^. 

If the condition of local thermodynamic equilibrium is 
satisfied, according to the definition, the auxiliary site 
attached to the system must be in equilibrium and then 
Eq.[Tl]must be satisfied. Conversely, if the linear relation 
in Eq. [Tl]is satisfied at the attached site, we would make 
sure that the condition of local equilibrium is fulfilled. 

In summary, we let a system evolve into a nonequilib- 
rium steady state. Then we couple an auxiliary site to a 
local region of the system. The coupling strength should 
be infinitesimal. We switch on the coupling and calculate 
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the occupation number at the site in steady hmit. We 
adjust the energy level at the auxiliary site to obtain 
a series of pairs (e^, rid). If we cannot fit Eq. [11] to the 
series of data points, we conclude that the local region 
is not in equilibrium. If we succeed, we then determine 
the local temperature and chemical potential from the 
fitting parameters. This is our method of using one aux- 
iliary site to calculate the local intensive thermodynamic 
parameters in a NESS. 

Next we give some comments on the possible applica- 
tions of it. Up to now solving the equation of motion 
is the only generally accepted way to reach a NESS in 
mesoscopic transport. In this approach we usually as- 
sume a central conductor connected to at least two in- 
finite reservoirs which are in equilibrium at initial time. 
The boundary conditions of the central conductor are 
given by the temperatures and chemical potentials of the 
reservoirs, which are assumed not to vary with time. The 
steady current through the central conductor can be cal- 
culated. The curve of current vs. voltage bias is plotted 
and compared with experimental results. However, some- 
times we are more interested in the distribution of local 
chemical potentials inside the central scattering region, 
which cannot be resolved from boundary conditions. Our 
method provides a way to calculate the local T and /i 
in the central conductor. The results can be compared 
with experiments. The local T and /i are more impor- 
tant than reservoir temperature and chemical potential 
in the theory of nonequilibrium thermodynamics, where 
the gradients of intensive variables are called thermody- 
namic forces and are the reasons of flows. 



IV. THE DISTRIBUTION OF CHEMICAL 
POTENTIALS INSIDE A DISORDERED CHAIN 

A chain connected to two reservoirs has been widely 
used to model the quantum point contact and the quan- 
tum wire fabricated in semiconductor heterostructures. 
If the two reservoirs are in different chemical potentials, 
a current will be driven through the wire. This is the 
Landauer-Biittiker picture for mesoscopic transport. In 
ballistic transport regime when the wire is clean, there 
is only contact resistance in the device. In other words, 
there are sharp voltage drops at the interface between 
the reservoirs and the wire, while the voltage is equal 
everywhere inside the wire. 

What is more interesting is a dirty wire, where we 
would expect the intrinsic resistance and a gradual drop 
of voltage through the wire. As in a macroscopic circuit, 
in a dirty wire the current can be treated as the result 
of the gradient of voltage inside the wire instead of the 
voltage bias of the reservoirs. The reservoirs are no more 
than the simulation of the voltage source in a circuit, 
while they are necessary for preparing the NESS. 

A disordered chain is obtained by proposing random 
onsite potentials e^, where j denotes the chain site. For 
simplicity we suppose that there is no correlation between 



ej at different sites. The average of ej is denoted as e. We 
suppose that Cj has a uniform distribution at the interval 
[e — where A denotes the disorder strength. 

Then the total Hamiltonian including the two reservoirs 
is written as 

kaa ka 
n n— 1 

+ '^A'y^i'^ + 9Y{c]aCj+i,a + h.c.) 

n 

+i7^c]^c,tc|cj^, (12) 

where Ckaa is the annihilation operator of the electron 
in the lead, a —tii denotes the spin and a — L,R de- 
notes the left and right lead. The leads are in thermal 
equilibrium at temperature and Tr, and in chemical 
potential fiL and fin respectively. The operator Cja- with 
j — 1,2- ■ - n denotes the annihilation operator at site j 
in the chain. The coupling strength between nearby sites 
in the chain is set to be g, and that between the leads 
and the chain is set to be gi. 

For general parameters, the above Hamiltonian is non- 
integrable due to the existence of Coulomb interaction. 
Several complicated numerical methods and approxima- 
tion schemes have been developed for solving this model. 
For simplicity, in this paper we will first solve the model 
without considering interaction by setting U = 0. Then 
we consider the effect of the interaction in self-consistent 
mean field approximation. Without interaction the elec- 
trons in different spin channels will transport indepen- 
dently and the spin index can be neglected. The simpli- 
fied Hamiltonian becomes 

H = ^kcljka + giY^icl^ci + cIj^c^ + h.c.) 

k.a—L.R k 
n n—1 

+ ^^3^% + 9Y{c]cj+i + h.c). (13) 

In the left and right electron reservoirs we employ the 
momentum basis and assume a constant density of states 
denoted by p. At the same time the reservoirs can also 
be simulated by the semi-infinite chains. In this sense 
the Hamiltonian of Eq. [13] describes an infinite chain 
with emphasized central sites, similar to the Hamilto- 
nian of Eq. [3] If we define F = pi^g^ as the energy unit 
(like what we do in single impurity Anderson model) and 
set g F, the coupling strength between central sites 
has the same amplitude as that between nearby sites in 
reservoirs. This corresponds to a comparatively strong 
coupling between central sites. In following text we al- 
ways set g = T. One could also take different values of g. 
But as g <^T the transmission spectrum has very sharp 
peaks, and one should be very careful when performing 
numerical integration routines. 



6 



A. The electron density and the current at site m 

The steady current and electron density at arbitrary 
site m in the disordered chain are got by the Keldysh 
technique, in which an adiabatic evohition is proposed 
with g and gi switched on little by little. In this model 
adiabatically switching on the coupling makes no differ- 
ence with a quench of coupling for the physical quantities 
in steady limit. In details we define the retarded Green 
functions as 

Gl^{t,t') = -t9it-t'){{c,{t),c]{t')}), (14) 

and the lesser Green functions as 

G<^{t,t')^^{c]ity4t)), (15) 

where i,j = l,2,---n denote the sites in the chain. The 
electron density and current at arbitrary site in the chain 
can be related to the simultaneous lesser Green functions. 

The freqency representation of Green functions is de- 
fined as 

G[;<(c.) = fdit- t')e-(*-*')G[:<(t,t'). (16) 



In frequency representations the retarded Green function 
satisfies the Dyson equation, i.e. 

^(a;) = G"''ioj) + G"'''ioj)j:'-{oj)G''iij), (17) 

where G^{uj) is a n x n matrix whose elements are G^ j (lu) . 
G^^{uj) is the free Green function matrix when g = gi = 
0, and its elements are 



1 



uj — ej + ii] 



(18) 



where 77 > is infinitesimal. S'"(a;) is the self-energy 
matrix, which comes from two sources. The first is the 
hopping energy between nearby sites in the chain, i.e. 



y — V — n 



(19) 



The second is the broadening of levels at the edge sites 
due to the coupling to leads, i.e. 



.9; 



2 \ " ^ 

bj - ek + iri 



(20) 



where T — pirgf. The other elements of the self-energy 
matrix are all zero. By sloving Eq. [17] we express the 
retarded Green function as 



^(w) = (G"'-(tj)-i - E'^(w)) 



(21) 



Here we need to calculate the inverse of a n x n matrix, 
which can be finished by computer. 



The lesser Green function in frequency representation 
is related to the retarded one in Keldysh formalism, i.e. 



,(22) 



where Gl<{uj) = 27 



1 



5[uj — efc) with a = 



L,R is the lesser Green function in left and right leads 
respectively. Here Pa = ^/{kBTa) and fXa are the inverse 
of temperature and the chemical potential in lead a re- 
spectively, and Gf j denotes the advanced Green function, 
which is the complex conjugation of the corresponding 
retarded Green function, i.e. Gfj 
uates 



(G^J*. Eq.[22]eval- 



where /q(w) 



Gf^j — 2ir {gi iGi j/l + Glj^G'^ jfji) , (23) 
1 



is the Fermi function in 



lead a. 

At last we perform inverse Fourier transformation and 
get the lesser Green function in time representation, i.e.. 



G< (i,0 



du)- 



-iu(t-t') 



27r 



-G<(c.). (24) 



The electron density at site m is defined as = [c^Cm), 
which can be expressed as Um = —iG"^ rn(0, 0) according 
to the definition of the lesser Green function. Substitut- 
ing Eq. [24] in, we have 



= ^ / dwG<_„(a;) 



(25) 



The current from site m to (m + 1) is expressed as 
Im^m+i = |(*(cLcm+i ^ h.c.)) , which Satisfies the 

electron number conservation law — — = Im-i^m ~ 

dt 

Im^m+i- Similarly it is related to the lesser Green func- 
tion as 



-25, 



Re / dujG. 



< 



(26) 



Eq. [2TI [23] [25] and [26] together give the complete proce- 
dure for calculating the electron density and current at 
arbitrary site. 

In numerical approach, we first use a random number 
generator to generate a group of onsite potentials ej . One 
should notice that the distribution of ej is uniform at 
the interval [e — tI- Then the electron density 

and current are calculated subject to this group of £j. 
Obviously repeating the calculation will not produce the 
same result, since ej is randomly generated. We are more 
interested in the average value of current, which is got by 
repeating the calculation for many times. In each time a 
new group of ej is generated. We repeat this procedure 
until the averaged current is convergent. 
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B. The onsite temperature and chemical potential 

We calculate the temperature and chemical potential 
at arbitrary site in the chain. This is done by attaching 
an auxiliary site to the site we want to measure. For 
example, if we want to measure the temperature T,„ and 
the chemical potential fj,m at site m, we will modify the 
Hamiltonian of Eq. [13] by including the extra terms 

Hmea = 9 {d) Cm + h.C.) + €dd) d, (27) 

where c„i and d denote the annihilation operators at site 
m and the auxiliary site respectively. Here ed denotes the 
energy level at the auxiliary site and is adjustable, and g' 
denotes the auxiliary coupling strength which should be 
infinitesimal. The whole system consisting of the chain 
coupled to two leads and the auxiliary site is described 
by the Hamiltonian H + Hmea, where H is defined in 
Eq.IH 

Again the Keldysh technique is employed to calculate 
the occupation number at the auxiliary site Ud = (d^d) 
corresponding to different ed- This should be done by 
adiabatically switching on g' after the chain described by 
H has been in nonequilibrium steady state. Due to lack 
of initial correlation, this is equivalent to switching on gi , 
g and g' simultaneously when we are only interested in 
the steady limit of the local quantity n^. 

Now the chain and the auxiliary site together are 
treated as the central conductor, which totally contains 
(n+ 1) sites. We define the retarded Green function con- 
cerning the auxiliary site as 

Gr,„+i(i,t') = -ie{t-t'){{c,{t)J\t')}), 

G;+i,,(t,t') = -ie{t-t'){{d{t),c]{t')}), 
GUi,n+iit,t') = -ieit-t'){{dit)JHt')}), (28) 

where i, j — 1,2, ■ ■ ■ n. Together with the Green functions 
defined in Eq.[T3]we get a (n+l) x (n+l) matrix. Similar 
to Eq. [21] it can be expressed as 

G'^(c.)= (g°'-(c.)-1-S'-(c.))"\ (29) 

where G'^^{lo) is the {n+l) x (n+l) diagonal matrix. The 
elements of G'^^ are defined in Eq. |TH] except for the one 

concerning the auxiliary site G^']-, . 

Now the self-energy matrix comes from three sources: 
the hopping energy between nearby sites in the chain; 
the coupling energy between the edge sites and the leads; 
and the coupling energy between m site and the auxiliary 
site. The corresponding elements are 

^^.+1 = ^^+1.. - 9 (30) 
with j — 1,2, ■ ■ ■ ,n ~ 1, 

Kl = Kn = -*r, (31) 



and 

The expression of the lesser Green function in Eq. [23] 
keeps invariant because the auxiliary site is not directly 
coupled to the leads. Especially we have 

G'n+i,n+i = 2ir (|G^j+i_iP/l + • (33) 

The occupation number at the auxiliary site is expressed 
as 

"rf = ^y" c'wG<+i_„+i(a;). (34) 

We have known that g' — > is the necessary condi- 
tion of local equilibrium at the auxiliary site. In nu- 
merical approach, the weak coupling limit is realized by 
self-adaptation algorithm. We begin with a finite g' , and 
reduce it by half in each loop until the final result Ud 
converges to a desired precision. 

At beginning of the numerical approach, we generate a 
group of onsite potentials e^. Then we in turn calculate 
the current through the chain and the occupation number 
at the auxiliary site with respect to different ed and to. 
When calculating the current, we turn off the coupling 
to the auxiliary site. One should notice that the destiny 
of an auxiliary site is to measure the temperature and 
the chemical potential at sites in the chain. One does 
not measure different sites simultaneously but rather one 
after another. We finally get a series of functions nd{ed) 
as m varying from 1 to n. We then calculate the average 
of temperature and chemical potential at site m with 
respect to different configurations of disorder. 



C. Transition from delocalization to localization as 
the disorder strength A increasing 

What we are studying is the transport through a meso- 
scopic chain. At particle-hole symmetry we set the av- 
erage of the Fermi energies in left and right leads to be 
equal to the average of the onsite potentials in the chain, 
i.e. (^L -I- ^i?)/2 = e = 0. The linear conductance of the 
chain in the absence of disorder (A — 0) reaches the uni- 

tary limit — . In the absence of disorder the transport 
h 

through the chain is well understood under the frame- 
work of ballistic transport. 

In the presence of disorder, however, we observe the 
gradual transition from delocalization to localization as 
the disorder strength increasing. This transition should 
be attributed to the Anderson localization. While the 
transition in our model is a "smooth" one, unlike the 
phase transition in thermodynamic limit, since the chain 
consists of only a few of sites. Strictly speaking, there is 
no disorder-induced phase transition in an infinite one- 
dimensional system, where an infinitesimal concentration 
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FIG. 3: The transmission spectrum —2gI{.eG^ at dif- 

ferent disorder strength. In this figure the chain length is set 
to be 50. We set V ^ - fJ-R ^ O.SF and Tl=Tr^ Q.OIF. 



of impurities will cause the exponential decay in the ex- 
tension of the wave function. In a chain with a few 
of sites, we see the gradual transition from a continu- 
ous band to sharp peaks and wide forbidden regimes be- 
tween peaks in the transmission spectrum as the disorder 
strength A increasing (see Fig. 3). This is understood as 
the transition from ballistic transport in the clean limit 
to the transport shuttled by the localized levels weakly 
coupled to the reservoirs. 

The feature of Anderson localization is clearly observed 
in the figure of averaged current (see the bottom figure 
of Fig. 4), in which the current shows an exponential 
decay as the length of the chain increasing. When we 
calculate the current, the onsite potentials are randomly 
generated according to I and A. For each generation 
of onsite potentials we obtain a different current. The 
averaged current is expressed as 



M 



(35) 



and the standard deviation of current is expressed as 
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(36) 



where Ij is the result of current corresponding to the j-th 
generation of onsite potentials and M is the total number 
of generations. 

The top figure of Fig. 4 shows the average and the stan- 
dard deviation of current varying with disorder strength. 
In clean limit as A = 0, the current approximately ap- 

2 

proaches to j^V and the standard deviation is exactly 
zero. As the disorder strength A increasing, the aver- 
aged current drops monotonously, while the standard de- 
viation first increases until its maximum and then drops 
towards zero. At the intermediate disorder strength, ct/ 
reaches its maximum. The strong fluctuation of current 
in this regime indicates the transition from delocalization 
to localization. 



D. Local equilibrium condition and the onsite 
temperature and chemical potential 



FIG. 4: The top figure shows the averaged current I and 
the standard deviation of current ai as a function of dis- 
order strength A. The voltage bias between two leads is 
V = HL — yiR = 0.05r, and the temperature in leads is 
Tl ^ Tr = O.ir. The length of the chain is set to be 10. 
We observe the localization-delocalization transition in the 
regime of intermediate disorder strength, where the deviation 
of current reaches its maximum. A large value of ai indicates 
a strong fluctuation of current, which often happens near the 
critical point of a phase transition. The bottom figure shows 
the averaged current as a function of the chain length. The 
voltage bias is = 0.2r. 



Local thermodynamic equilibrium is an important con- 
cept in the study of nonequilibrium thermodynamics. 
Generally the local equilibrium condition is only fulfilled 
in the systems not far from equilibrium. In the case of a 
chain connected to two leads, it requires that the temper- 
ature difference and the voltage bias of the two leads are 
small. We have explained how to use an auxiliary site to 
check if the local equilibrium condition at arbitrary site in 
the chain is fulfilled by calculating the occupation num- 
ber at the auxiliary site as a function of the energy level. 
In Fig. 5 we plot ln(l/n(j — 1) as a function of under 
different voltage bias, which is defined as V = — fJ-R- 
According to Eq. 1111 the local equilibrium condition is 
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FIG. 5: An auxiliary site is attached to the third site in a 

chain of length 10 coupled to two loads. Whether the con- 
dition of local equilibrium is fulfilled at the third site is de- 
termined by the relation between the occupation number 
and the energy level td at the auxiliary site. This figure 
shows ln(l/nd — 1) as a function of ed at different voltage bias 
1/ = 0, 0.05r, 0.2r, O.SSr and O.Sr. Here we set e = and the 
disorder strength A = 0. The temperature in left and right 
leads is Tl = = O.IF. One should notice that the curves 
corresponding to F = and V — 0.05T have strong overlap 
with each other and are difficult to be distingished. Obviously 
when the whole system is in thermodynamic equilibrium as 

V = 0, the site we are studying is exactly in local equilibrium, 
as shown in the straight line titled = 0. The slope of this 
line is easily found to be lO/T, which is just the inverse of local 
temperature. The line intercepts the axes at original point, 
indicating the local chemical potential is zero. If a voltage 
bias as small asV = COST is applied to the system, the func- 
tion ln(l/nd — 1) keeps approximately a linear function and 
then the condition of local equilibrium is still fulfilled. But 
as voltage bias increasing, the function ln(l/nd — 1) obviously 
deviates from a linear function, as shown in the curve titled 

V = O.SF. 
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FIG. 6: The value of ln(l/n£j — 1) is plotted with respect to ed 
at different sites in a disordered chain of length 10. The left 
and right lead temperature is Tl ~ Ta = O.IF. The voltage 
bias is set to be F = 0.05F. The disorder strength is A = F. 



FIG. 7: The distribution of onsite chemical potentials inside 

a chain of length 10 at different voltage bias V and disorder 
strength A. Strong fluctuation of onsite chemical potentials 
is observed at finite voltage bias in the presence of disorder. 




The site m 

FIG. 8: The distribution of averaged onsite chemical poten- 
tials. Here we set the voltage bias to be F = 0.05F and 
the disorder strength to be A = F. By taking average with 
respect to configurations of disorder we observe the regular 
distribution replacing strong fluctuations. 



fulfilled if this function is a straight Une. As shown in 
Fig. 5, when V is exactly zero or very small the hmction 
of ln(l/nd — 1) is in fact linear. While as V increasing 
to as large as O.SF the function obviously deviates from 
a linear function. Then the function nd{&d) cannot be 
regarded as a Fermi function at large voltage bias, and 
the concepts of temperature and chemical potential are 
nonsense. 

An interesting question arises as to whether the local 
equilibrium condition is less stringent if the auxiliary site 
is simultaneously coupled to a large number of sites in 
the chain (and not only to a given site). We study a 
model in which the auxiliary site is coupled to up to 
100 sites in a chain of length 200. Under finite voltage 

bias, the deviation of the function ln( — - — - — 1) from a 

nd{td) 

straight line is found to be indepedent to the number of 
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sites measured simultaneously, in opposite to the results 
observed in an equilibrium system [25j. 

At small voltage bias, the function \ii{l/nd — 1) is lin- 
ear at arbitrary site inside the chain. The slope of the 
function is explained as the inverse of local temperature, 
and its value at = as (— where ^ is the local 
chemical potential. Then we can extract the local tem- 
perature Tjyi — 1/ Pm and the local chemical potential /i™ 
at each site m. One may have doubt on the validation of 
this method. Since what we are interested in is the distri- 
bution of onsite temperatures and chemical potentials in 
nonequilibrium, when they are different from the reser- 
voir ones. However at finite voltage bias, the condition of 
local equilibrium is only fulfilled approximately. Strictly 
speaking, the function \n{l/nd — 1) deviates from a linear 
function a little bit even at small voltage bias. Then how 
do we know that the discrepancy between onsite chemi- 
cal potentials and reservoir ones if there is any is not the 
result of the deviation of ln(l/n£; — 1) from a linear func- 
tion? This suspicion is canceled by plotting the function 
ln(l/nd — 1) at different sites simultaneously (see Fig. 6). 
In Fig. 6, we see that subject to a finite voltage bias and 
disorder strength the functions of ln(l/nd — 1) at differ- 
ent sites are parallel to each other. They have obviously 
different intercepts with the axes and at the same time 
the same slope which is just the inverse of the reservoir 
temperature. This is a strong evidence that the deviation 
of ln(l/nd — 1) from a linear function is very small com- 
pared to the difference between onsite potentials. Then 
our method resolves the onsite chemical potentials to a 
good extent. 

We study the distribution of onsite chemical potentials 
at different voltage bias and disorder strength. Even if 
the result depends on the configuration of onsite poten- 
tials and is not repeatable when A 7^ 0, it shows some 
features which can be attributed to V and A (see Fig. 7) . 
In equilibrium as = 0, the onsite chemical potential is 
zero throughout the chain whether there is disorder or 
not, coinciding with the fact that the chemical potential 
in an equilibrium system should be equal everywhere. 
When there is no disorder, the onsite chemical potential 
keeps zero everywhere even at finite voltage bias. Be- 
cause without disorder the resistance of the device comes 
only from contact resistance, and there is no drop of 
voltage inside the chain. In the presence of both volt- 
age bias and disorder, we see strong fluctuation of onsite 
chemical potentials. As disorder strength increasing, the 
fluctuation becomes stronger. This is contrary to the 
familiar phenomena in macroscopic circuits where we ex- 
pect a monotonic drop of voltage through the resistance 
wire. This can be explained by the mesoscopic nature 
of the chain, where the movement of electrons must be 
described by quantum mechanics. In localization regime 
with intermediate A, the position of the electron in the 
chain is constrained to a small region, then its kinetic en- 
ergy has strong fluctuation due to uncertainty principle, 
leading to the strong fluctuation in the onsite chemical 
potentials. 



The abnormal fluctuation of onsite chemical potentials 
is then canceled by taking the average with respect to 
configurations of disorder. Similar to the averaged cur- 
rent in Eq. 1351 we define the averaged chemical potential 
as 



M 



M 



(37) 



where /i^ denotes the chemical potential at site m in the 
j-th. generation. In a real mesoscopic wire, the configura- 
tion of disorder is generally not controllable and unknown 
to the simulators. By taking average with respect to the 
disorder configurations, one could obtain some universal 
feature about the distribution of onsite chemical poten- 
tials. This is a well-established method in studying disor- 
dered systems. As shown in Fig. 8, the averaged chemical 
potential drops monotonously through the chain, recov- 
ering the feature in macroscopic resistance wires. 



E. The effect of interaction 

We have studied the distribution of onsite chemical po- 
tentials in a non-interacting disordered chain. A natural 
question arises as to what is the influence of electron- 
electron interaction to our results. Next we consider the 
complete Hamiltonian in Eq. [T^] consisting of the onsite 
Coulomb interaction. We combine it with the measuring 
term 



Hmea = ff' ^(J^Cmo- + h.C.) + ^ d^da 



(38) 



The two spin channels are considered at the auxiliary site, 
while the onsite interaction is not. Because if we consider 
the Coulomb energy, there would be four states with dif- 
ferent energies at the auxiliary site. This increases the 
difficulty in deciding the local temperature and chemical 
potential. 

Again the Keldysh techniques are employed to calcu- 
late the occupation number at the auxiliary site as a func- 
tion of Ed in NESS. The Green function can be expressed 
in a perturbative expansion of U. The Green functions 
of zeroth order in U have been given in above calcula- 
tions, i.e., in Eqs.l^andB^ We take the self-consistent 
mean field approximation by summing up all the dia- 
grams in the perturbative expansion in which Hartree 
type self-energies are inserted. This is finished by replac- 
ing the onsite potentials ej by ej + Urij in the expression 
of G°''(a;) in Eq. [THl hence in expressions of G'^{uj) and 
G^(aj) in Eqs. [H] and [231 Here Uj is the occupation 
number at site j and satisfies the self-consistent equation 



1 

27ri 



dujGf.j{uj). 



(39) 



Above self-consistent equations can be solved in the it- 
erative method as the interaction strength U is small. 
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FIG. 9: The distribution of onsite chemical potentials at dif- 
ferent interaction strength U and disorder strength A. The 
chain is set to be at particle-hole symmetry with e = —U/2. 
The voltage bias is V = O.OSF, and the left and right lead 
temperature are both O.IT. Without disorder the interac- 
tion has no effect on the onsite chemical potentials (see the 
two curves titled A = 0). In the presence of disorder, the 
electron-electron interaction will suppress the fluctuation of 
onsite chemical potentials caused by localization of electrons. 

We plot the distribution of onsite chemical potentials 
at the particle-hole symmetric point where e = —U/2 (see 
Fig. 9). We find that the interaction will suppress the 
fluctuation of onsite chemical potentials in the presence 
of disorder. This is similar to the phenomenon that in a 
closed system the interaction will cause equilibration of 
the system from arbitrary initial states. 
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